Joint modeling of longitudinal measures of pneumonia and time to convalescence among pneumonia patients: a comparison of separate and joint models

Background Globally, pneumonia is the leading cause of children under age five morbidity and mortality with 98% of deaths in developing countries. Objective This study aimed to identify the determinants of longitudinal measures of pneumonia and time to convalescence or recovery of under five admitted pneumonia patients at Felege Hiwot Referral Hospital, Bahir Dar, Ethiopia. Methods A prospective cohort study was conducted among a randomly selected sample of 101 pneumonia patients using simple random sampling who were on follow up from December 2019 to February 2020. A Linear mixed effect model were used for the longitudinal outcomes and joint model for modeling both longitudinal and time to event outcomes jointly respectively. Results The significant values of shared parameters in the survival sub model shows that the use of joint modeling of multivariate longitudinal outcomes with the time to event outcome is the best model compared to separate models. The estimated values of the association parameters: − 0.297(p-value = 0.0021), − 0.121) (p-value = < 0.001) and 0.5452 (p-value = 0.006) indicates association of respiratory rate, pulse rate and oxygen saturation respectively with time to recovery. The significant values show that there is an evidence to say that there is a negative relationship between longitudinal measures of respiratory rate and pulse rate with time to recovery and there is positive relationship between longitudinal measures of oxygen saturation with time to recovery. Variables age, birth order, dangerous signs, severity and visit time were significant factors on the longitudinal measure of pulse rate. The significant factors related to longitudinal measures of oxygen saturation were birth order, severity and visit. From this we can conclude that birth order, severity and visit were significant variables that simultaneously affect the longitudinal measures of respiratory rate, pulse rate and oxygen saturation of patients at 5% level of significance. Conclusion Results of multivariate joint analysis shows that severity was significant variable that jointly affects the three longitudinal measures and time to recovery of pneumonia patients and we can conclude that patients with severe pneumonia have high values of respiratory rate and pulse rate as well as less amount of oxygen saturation and they need longer time to recover from the disease.


Background
Pneumonia is described as the inflammation of parenchymal structures of the alveoli and the bronchioles (lungs) [1]. Community-acquired pneumonia (CAP) is an infection that begins outside the hospital and/or diagnosed within 48 hours after admission to the hospital. Whereas, hospital-acquired pneumonia occurs in more than 48 hours after admission and without any antecedent signs of infection at the time of hospital admission [2,3].
Pneumonia is usually caused by infection with bacteria or viruses and bacteria are most common causes of (CAP) with Streococcus pneumonia isolated about 50% of cases. In children about 15% pneumonia cases a number of drug-resistant versions of the infections are more common, including drug resistant Streptococcus pneumonia and Methicillin-resistant Stapylococcus aurous [4,5]. The burden of medical response to pneumonia has significant challenges. Besides drug resistance to the bacteria, comorbid conditions like Malaria, TB, Sickle cell anemia, HIV/AIDS and risk factors like lack of exclusive breast feeding, alcoholism, smoking etc. commonly appear in pneumonia patients which leads to define the severity and risk scores of the disease in which used for clinicians to make care self-site decision as in-patients or out-patients [6].
Estimates from the WHO suggest that pneumonia is responsible for 20% of deaths in the under-five age group, leading to 3 million deaths per year [7]. In Africa especially in sub-Saharan Africa, by 2013 pneumonia was the second leading cause of child mortality that accounts a million child death about 15.8% of total deaths in the region [8]. A report by UNECIF indicated that 132,000 under five children killed by pneumonia in Congo which is the second cause of death next to malaria [9,10] and Kenya accounted the highest number of under-five mortality due to pneumonia which accounts about 16% of total deaths among 15 East Africa countries [11]. Pneumonia kills up to 5 million children under the age of 5 years annually in developing countries [12].
In Ethiopia, pneumonia is a leading single disease killing under five children and it contributes about 18% of all cases (370,000) of under five deaths compared to diseases like diarrhea, AIDS, malaria and measles every year [13,14]. Under five pneumonia is commonly measured through physiologic parameters (temperature, pulse rate, blood pressure, respiratory rate, and oxygen saturation) and the performance of (TCS) is decided through the longitudinal measures of those parameters [15].
Several cross-sectional studies have used scoring systems to summarize the level of symptoms within a cohort at fixed time-points following CAP [16,17]. However the understanding of which predictor affects length of hospital stay has been hampered by lack of longitudinal studies. Recent studies provide insight on the background and clinical predictors of mortality and survival of pneumonia patients among children aged under 5 years [18]. These studies did not consider the true and unobserved effects of longitudinal measures of physiologic parameters which correlates with recovery time to determine the length of hospital stay for under five pneumonia patients. In this study, joint model of multivariate linear mixed model and cox PH model was used to find significant factors of longitudinal measures of pneumonia (respiratory rate, pulse rate and oxygen saturation) and time to convalescence jointly.

Study area and period
The data for this study was collected from FHRH, Bahir Dar, Ethiopia from 12th December 2019 to 30th February 2020. Bahir Dar is the capital city of Amhara National Regional state. It is found in north western Ethiopia and is 565 km from Addis Ababa. This hospital serves as referral hospital for the people who came from different surrounding areas.

Study design and sampling
A prospective cohort study design was conducted from 12th December 2019 to 30th February 2020. The study population was all selected under five pneumonia patients who were on treatment and follow up at FHRH from 12th December 2019 to 30th February 2020 and who full fill the inclusion criteria. The samples were the number of children bounded in the inclusion criteria with in the study period since the study is a case study carried out within 3 months follow up period. Therefore, all the confirmed cases of severe pneumonia described as per WHO criteria, aged from 2 to 59 months were included in the study. Hence, the final sample size estimated for this study was 101 patients. To select the study participants, simple random sampling method using table of random numbers was employed.
severe pneumonia have high values of respiratory rate and pulse rate as well as less amount of oxygen saturation and they need longer time to recover from the disease.

Variables in the study
Three longitudinal outcomes (Respiratory Rate in bpm, Pulse Rate in bpm and Oxygen saturation in mm Hg) and a survival outcome (time to convalescence or recovery in hour) were considered as dependent variables in this study.

Data collection procedure
The longitudinal and the survival data containing the socio demographic and home based information were collected using primary data collection method by face to face interview of their care givers using well-structured questionnaires. In addition, the data containing clinical information found from their charts were considered. Both primary and secondary data were collected by trained pediatrician and statistician.

Eligibility criteria
The inclusion criteria was children 2-59 months of age with their care givers (mothers or not) and admitted at pediatric ward by community acquired pneumonia during the study period. Whereas the exclusion criterion was children admitted at the hospital by disease other than pneumonia, pneumonia patients below 2 months and above 59 months and patients with incomplete medical records.

Data analysis
In this study, a longitudinal data on the three measures of pneumonia (RR, PR and OS), recovery time of underfive pneumonia patients for the survival data, and sociodemographic factors, home based factors, child nutritional status and child illnesses at the base line were considered. The data were coded, entered and edited using SPSS version 26 and the analysis was done using SAS 9.4 and R software and the statistical decision was made at 5% level of significance.

Survival data analysis
Survival data analysis is a class of statistical method which used to analyze data in which the time(usually measured in days, weeks, months or years) until the event (usually death, disease incidence, relapse from remission, recovery) is of interest [19]. Cox proportional hazards model of the survival analysis was used to estimate the length of time to recover from pneumonia and to identify factors related to time to recovery [20]. The Cox model is defined as:

Longitudinal data analysis
A longitudinal study is statistical analysis of an observational research method in which response variable is measured repeatedly over time and those measurements taken from the same subject are correlated [21]. Longitudinal response may arise when measurements taken on the same subject or when measurements taken on related subjects. In both cases, the responses are likely to be correlated [22].

Linear mixed effects model
The random effects contains subject specific random effect and are directly used in modeling the random variation in the dependent variable at different levels of the data. Before considering the multivariate linear mixed model, it is better to identify the covariates which have significant effect on the mean change of RR, PR and oxygen saturation measurements over time using LMM [21].
Let y ijk represent the j th observation of the k th outcome variable for the i th subject, where: The vector (y 1ik + y 2ik + . ……. . + y nik ) T represents the n ik observations of the k response variable from the i th subject and vector (y 1k , y 2k , …. . y nik) T represent the N k observation for the k th response variable across all response variables and subjects, finally the vector (y 1 , y 2 , y 3 , ………y n ) T represents the observations across all response variables and subjects. In the context of modeling the response variables, the linear mixed effect model for each response variable of subject i, taken at time t, can be specified by [23].
where μ k (t) is the average evolution of the k th response over time and it is a function of fixed effect. The subject specific random intercepts α ik and slopes b ik (t) describe how the subject specific profiles deviate from the average profile for the k th response.

Joint modeling of multivariate longitudinal with time to event outcome
In this study three correlated and longitudinally measured response variables were considered which can be jointly modeled with time to event outcome. The separate and the joint models assume that the longitudinal sub model has the form similar to the conventional i = 1, 2, … … 101, j = 1, 2, … … .n i and also, K = 1, 2, 3.
(2) linear mixed effects model while the survival model in the joint model includes a latent association function w i (t) [24]. Maximum likelihood approach was used to estimate the parameters for both longitudinal and survival sub models.

Ethical consideration
This study was carried out in the location where the approval was obtained from the ethical review committee of College of Health Sciences, Bahir Dar University, and permission for data collection was obtained from Felege Hiwot Specialized Referral Hospital Management. There were no risks due to participation in this research project, and the collected data were used only for this research purpose. The study compiled with the principles set forth in the Declaration of Helsinki (1964) and all of its subsequent amendments. The written informed consent was obtained for caregivers of each patient prior to the data collection and all information collected from each caregivers was treated with complete confidentiality.

Results
The study revealed that, the median recovery time of pneumonia patients admitted at FHRH was 72 hours with minimum and maximum recovery time of 18 hours and 96 hours respectively. Out of the total sampled pneumonia patients, 90 (89.1%) were recovered from pneumonia. When we fit the cox proportional hazards model using the candidate variables: residence, birth order, age of mothers, education of mothers, danger signs, cooking place, comorbidity and severity were significant factors affecting time to recovery of pneumonia patients at 5% level of significance (Table 1).

Separate analysis of longitudinal data
In this study, three longitudinally measured response variables of pneumonia patients were considered. The linear mixed model was used for all the variables; pulse rate, respiratory rate and oxygen saturation of patients. The study was started by exploring the mean and variance structure of those longitudinally measured response variables. The three longitudinal measures of pneumonia were approximately measured every 6 h a day from admission up to hospital discharge of under-five pneumonia patients. All of 101 sampled under five admitted pneumonia patients were at risk of pneumonia up to the third visit time (t = 12 hour), this tells that, for this study, the minimum follow-up time at which the patient get the event of recovery was the third visit (t = 18 hours) and the number of patients getting the event increases, whereas the number of patients with at risk of pneumonia decreases through visit time.
The study also revealed that, the average values of RR and PR decrease, whereas the average values of oxygen saturation increase through the visit time. At the end of the follow up, the overall average values of RR, PR and Oxygen saturation were 50.55 bpm, 131.20 bpm and 90.18 mmHg with standard deviation of 12.55 bpm, 27.37 bpm and 6.11 mmHg respectively (Table 2).
Checking assumptions of the data is the first step in analyzing longitudinal data. Normal QQ plots in Fig. 1 shows that, the data for the three longitudinal outcomes were approximately normally distributed and then it is better to proceed to the next steps of the analysis.

Multivariate analysis of longitudinal data
MLMM was fitted using three longitudinal measures of pneumonia (RR, PR and oxygen saturation) for underfive admitted pneumonia patients (Table 3). At 5% level of significance; marital status of mothers, smoking exposure of patients, breast feeding, severity, cooking place, comorbidity and visit time were significant factors related with longitudinal measures of RR. Age, residence, birth order, comorbidity, danger signs, vaccination, severity and visit time were significantly related with longitudinal measures of pulse rate. The variables that significantly related with longitudinal measures of oxygen saturation were; age at the base line, residence, comorbidity, danger signs, age of mothers, severity and visit time. Variables; severity, visit time and comorbidity were simultaneously  The random part of MLMM shows the variance and covariance between rate of change and baseline values for the three longitudinal measures of pneumonia (RR, PR and oxygen saturation) were significantly different from zero which tells the existence of a relationship between a patients baseline standing between outcomes, rate of change between outcomes as well as, between baseline standing of one outcome and rate of change of the other outcome through follow-up time.

Joint modeling of multivariate longitudinal data and survival data
In the previous sections; determinants of the multivariate longitudinal measures of pneumonia as well as determinants of time to recovery of under-five admitted pneumonia patients were identified. The results of joint model analysis for multivariate longitudinal and survival data found in the Table 4, contains multivariate longitudinal and survival sub models. In the random part of MLMM, estimates of variance and covariance were different from zero, shows the existence of correlation between intercepts of outcomes, between rate of changes of outcomes and correlation between rate of change and baseline values of the three longitudinal measures of pneumonia (RR, PR and oxygen saturation). Based results of Table 4, the average RR, PR and oxygen saturation of under-five pneumonia patients admitted at FHRH were 47.2660 bpm, 146.7431 bpm and 87.29 mmHg respectively when all categories are at their reference group. As age of patients increased by 1 month, the average RR and PR were significantly decreased by 0.38 bpm and 1.01 bpm respectively. Whereas, age was not a predictor of oxygen saturation. Coming from urban residence increases the average RR and PR by 1.70 bpm and 1.26 bpm respectively, whereas it lowers the average oxygen saturation by1.01 mmHg as compared with rural residency, keeping other variables constant. Being first child significantly rises the average RR and PR by 1.59 bpm and 3.09 bpm respectively; whereas it lowers the average oxygen saturation by 1.72 mmHg as compared with being second or above child; other variables held constant. Being non-exposed by smoking lowers the average RR by 2.27 bpm as compared with patients exposed by smoking; keeping other variables constant, but had no information about PR and oxygen saturation.
Being non-comorbid significantly lowers the average RR and PR by 3.98 bpm and 3.64 bpm respectively, while it rises the average oxygen saturation by 2.33 mmHg as compared with being comorbid, keeping other variables constant. Having sever pneumonia at the baseline increases the average values of RR and PR by 5.46 bpm and 1.30 bpm respectively, whereas it lowers the average oxygen saturation by 1.03 mmHg as compared with those having non-sever pneumonia, other factors held constant. Having literate mother increases the average oxygen saturation by 2.70 mmHg as compared with those from illetrate mothers, held other variables as constant. Cooking food inside the living room lowers the average oxygen saturation by 2.11 mmHg as compared with those whose parents cook their food out of living room, keeping remaining factors constant. A unit increase in the number of visits lowers the average RR and PR by 0.19 bpm and 0.16 bpm respectively, whereas it rises the average oxygen saturation by 0.90 mmHg, keeping other predictors constant.
Getting vaccination lowers the average PR by 8.59 bpm as compared with unvaccinated by remaining other variables constant. Feeding exclusive breast within first 6 months decreases the average RR by 1.85 bpm as compared with no breast feeding. The estimated hazard ratio of patients from urban area relative to patients from rural area was 0.61 indicates, patients from urban residence were 0.547 times less likely to recover from pneumonia than patients from rural residence, other variables held constant. Patients without comorbidity were about 2.296 times more likely to experience the event of recovery compared to patients without comorbidity. Patients at the first birth were 0.284 times less likely to get the chance of recovery compared to patients at the second and above births, keeping other variables constant. As age of mothers increase by 1 year, experiencing the event of recovery increases about 2.462 times, other variables held constant. Exclusively breast feed patients within first 6 months of life were about 4.06 times more likely to get recovery as compared with patients having no breast feed.
Patients with severe pneumonia were about 0.206 times less likely to experience the event of recovery compared to patients with non-sever pneumonia, keeping other variables constant. The estimated values of association parameters γ_1 = − 0.297 (p-value = 0.0021), γ_2 = − 0.121 (p-value< 0.001) and γ_3 = 0.545 (p-value = 0.006) indicates; RR and PR were negatively associated with time to recovery, whereas oxygen saturation was positively associated with time to recovery of under-five admitted pneumonia patients.
Model comparison: The multivariate longitudinal submodel was consistent with the results from the multivariate longitudinal analysis of RR, PR and oxygen saturation. The differences in magnitudes of the parameter estimates were negligible and there were some parameter difference in terms of statistical significance in separate MV longitudinal and separate survival model. But, longitudinal sub-model had narrow confidence interval which indicates that standard error is small for all significant predictors as compared to separate model in MV longitudinal and survival model. When evaluating the overall performance of both the separate and joint models in terms of model parsimonious and goodness of fit, the joint model was preferred as it has smaller standard error than the separate model. This result also supports the study done by [25,26].
As Table 4 revealed, under MV joint model, estimate of the association parameters in the survival sub model was significantly different from zero (γ_1 = − 0.297, γ_2 = − 0.121 and γ_3 = 0.5452), this indicates that three longitudinal outcomes were correlated with time to recovery of under-five admitted pneumonia patients supported by [27][28][29], stats that the longitudinal and survival data are correlated. The joint model was more parsimonious fit than the separate model. Therefore, the joint model found preferable and parsimonious to fit the data  better than the separate one [24] when the association parameter of the joint model is significant. Therefore, the final model for this study was joint model of MLMM and cox PH model.

Discussion
The general objective of this study was identifying the determinant factors jointly affecting longitudinal measures of pneumonia (RR, PR and oxygen saturation) and time to recovery of under-five admitted pneumonia patients at FHRH, Bahir Dar, Ethiopia and the discussion was made based on results of Table 4.
The result reveals that about 89.1% of under-five patients were recovered from pneumonia with a median time of 72 hours (3 days) which took shorter recovery time compared to results of the study done by [2,16,30], whereas it is longer recovery time compared to results of the study done by [18,31,32]. The difference can be due to the difference in explanatory that we used, type of hospital etc. Age has significant effect on the two longitudinal measures of pneumonia (RR and PR), but had no information about longitudinal measures of oxygen saturation. When age of patients increase, the RR and PR measures decrease for under-five admitted pneumonia patients. This indicates that, lower level of pneumonia are found for increased age of patients. This was in line with results of the study conducted using nonlinear mixed model by [16]. Unlikely, using binomial logistic regression [33] found that age had no significant effect on measures of pneumonia. This requires further investigation to reach a decision in the effects of age on CAP.
Urban residency significantly increases the average values of RR and PR, whereas it decreases the average values of oxygen saturation of under-five admitted pneumonia patients. This indicates that, urban residency was significantly associated with the risk of pneumonia. This contradicts with results of the study done by [33] using binomial logistic regression. Based on the two contradicted ideas, we can suggest that in our country most of the people living in urban area have not their own living house and they live within a crowded room by using as living room and cooking room which is difficult to treat children and to gain fresh air. This makes children to be highly vulnerable to pneumonia compared to patients from rural residence.
Smoking exposure increases the average values of RR and PR again it lowers the average oxygen saturation. This shows, smoking exposure related with increased level of pneumonia. This also agreed with results of the study conducted by [33] using binomial logistic regression for longitudinal data.
Having literate mothers increases the average oxygen saturation of under-five admitted pneumonia patients. This consides with results of the study conducted by [16]. Cooking food inside the living room lowers the average values of oxygen saturation which relates with high risk of pneumonia. As the follow-up time goes, the average values of RR and PR decreases, whereas the average values of oxygen saturation increases through visit time for under-five admitted pneumonia patients, which indicates effectiveness of treatment to lower pneumonia. In the survival sub model; Variables of urban residence, feeding exclusive breast within 6 months, first birth, non-danger sign and severity were significantly associated with recovery time of under-five admitted pneumonia patients. This was consistent with results of the study conducted by [31]. Increasing age of mothers increases the chance of experiencing the event of recovery (p-value = 0.026). This consides with results of the study conducted by [32]. The difference in the degree of significance may come from the difference in the variables as well as the model we used.
Exclusive breast feeding with in the first 6 months of life increases child survival by reducing the length of hospital stay. This supports results of the study done by [32]. The association parameters were significant indicates the significance of relationship between longitudinal measures of pneumonia (RR, PR and oxygen saturation) and time to recovery of under-five admitted pneumonia patients. This is in line with results of [34,35]. Higher values of average RR and PR as well as lower values of average oxygen saturation were related with longer recovery time (high risk of pneumonia). This was consistent with results of the studies done by [18,32].

Conclusion
In this study, a joint model of multivariate longitudinal changes of respiratory rate, pulse rate and oxygen saturation with time to recovery of under-five admitted pneumonia patients was discussed. Out of the total sampled pneumonia patients 90 (89.1%) were recovered from pneumonia and the median recovery time was 72 hours. When evaluating the overall performance of both the separate (MLMM and cox PH model) and joint model in terms of model parsimonious, goodness of fit and the statistical significance of association parameters, the joint model performs better than the separate models. As a result, we concluded that the joint model was preferred for simultaneous analyses of repeated measurement and survival data. From results of the study, we can conclude that patients from urban area, borned at the first birth, having comorbid status, age of mother, exclusive breast feeding and having sever pneumonia have high levels of respiratory rate and pulse rate, whereas lower levels of oxygen saturation and which increases the risk of pneumonia. Patients with high levels of respiratory rate and pulse rate as well as low values of oxygen saturation requires longer recovery time of under-five admitted pneumonia patients. To improve child survival, the health professionals and community should be responsible for post ponding child birth and marriage.